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The  availability  of  real-time  information  on  sea  ice 
conditions  in  ice  covered  seas  has  always  been  impor¬ 
tant,  not  only  to  strategic  military  operations,  but  to  the 
economies  of  those  countries  that  border  the  Arctic  and 
its  marginal  seas.  Knowledge  of  the  thickness  and 
movement  of  sea  ice  as  well  as  the  locations  of  open 
water  is  required  for  traversing  the  Arctic  whether  in  a 
drill  ship,  in  a  cargo  vessel  or  in  an  ice  strengthened 
ship  such  as  a  Coast  Guard  ice  breaker. 

Forecasting  these  conditions  is  a  difficult  task  at 
best.  The  ice  and  snow  that  cover  the  cold  Arctic  Ocean 
are  highly  variable  on  short  time  scales,  such  as  days  to 
weeks,  and  longer  time  scales  of  years  to  decades.  This 
variability  in  the  sea  ice  cover  is 
due  to  a  combination  of  dynamic 
and  thermodynamic  effects. 

Surface  stresses  on  the  top  and 
bottom  of  the  ice  cause  the  move¬ 
ment  of  sea  ice,  or  ice  drift,  as  well 
as  the  deformation  of  the  ice. 

Ideating  and  cooling  from  the 
atmosphere  and  the  ocean  are 
largely  responsible  for  the  growth 
and  decay  of  sea  ice,  in  combina¬ 
tion  with  the  ice  motion. 

Sea  ice  has  a  strong  seasonal  variability.  The 
thinnest  sea  ice  and  largest  amount  of  open  water 
appear  in  the  summer  months  from  June  to  September. 
Ice  begins  to  grow  back  in  the  fall  and  builds  to  a  max¬ 
imum  thickness  in  the  late  winter  and  early  spring, 
March-April.  Many  of  the  marginal  seas,  such  as  the 
Barents  and  Greenland  Seas  are  nearly  ice  free  in  the 
summer. 

Other  marginal  seas  like  the  Bering  Sea  and  the  Sea 
of  Okhotsk  are  completely  ice  free  in  the  summer. 
These  seasonal  patterns  vary  inter-annually  as  well. 
This  variability  is  often  represented  by  a  seesaw  effect 
when  one  part  of  the  Arctic  experiences  a  "mild"  ice 
year  while  another  part  of  the  Arctic  has  an  increase  in 


ice  extent  and/or  ice  thickness.  A  general  trend  overly¬ 
ing  this  inter-annual  variability  has  been  seen  by  sever¬ 
al  scientists  who  have  examined  long  records  of  satel¬ 
lite  data.  The  trend  is  for  an  overall  decrease  in  ice 
extent  from  the  late  1970s  through  the  mid-1990s 
(Parkinson  et  al.,  1999;  Cavalieri  et  al.,  1997; 
Johannessen  et  al.,  1995).  Whether  or  not  this  decrease 
will  continue  and  what  effect  it  may  have  on  Arctic  eco¬ 
nomics,  operations  and  climatic  change  is  a  topic  of 
great  international  interest. 

Sea  ice  forecasts  usually  focus  on  short  time  scales 
of  5-7  days.  In  order  to  provide  an  accurate  prediction, 
forecast  systems  are  most  often  based  on  a  combination 
of  models  and  data.  Modeling  sea 
ice  can  be  a  difficult  problem,  as  it 
exists  in  many  different  forms 
(Figure  1). 

It  can  appear  as  a  field  of  dis¬ 
jointed  ice  floes,  as  a  level,  contin¬ 
uous  field  of  ice  or  as  a  landscape 
of  ice  hills  and  ridges.  Interfacial 
stresses,  from  the  atmosphere 
and  the  ocean,  interacting  with 
natural  coastal  boundaries  cause 
the  ice  to  be  in  a  nearly  constant  state  of  motion  and 
deformation  resulting  in  the  formation  of  ice  rubble, 
ridges,  leads  and  ice  floes.  Most  sea  ice  models  make 
the  assumption  that  ice  exists  as  a  "continuum".  The  ice 
is  modeled  as  having  some  combination  of  viscous, 
plastic,  or  elastic  material  characteristics.  The  way 
that  the  ice  behaves  as  a  material  when  acted  upon  by 
stress  is  called  ice  rheology.  Ice  within  a  model  grid 
cell  is  defined  as  having  a  certain  thickness  and  cover¬ 
age,  also  defined  as  a  percent  of  ice  concentration. 
Thickness  and  concentration  are  allowed  to  vary  in 
time  depending  upon  the  atmospheric  and  oceano¬ 
graphic  conditions. 

The  Arctic  presents  a  hostile  environment  for  those 
trying  to  gather  data  on  the  atmosphere,  ice  or  ocean. 


Modeling  sea  ice  can  be  a 
difficult  problem,  as  it  exists  in 
many  different  forms. 
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Figure  1.  Different  forms  of  sea  ice:  a)  icefloes,  b)  a  level  field  of  solid  sea  ice  penetrated  by  the  sail  of  the  submarine  USS 
Hawkbill,  SCICEX  ‘99  and  c)  ridged  ice.  Figures  are  courtesy  of  a)  the  Canadian  Ice  Services,  b)  the  Arctic  Submarine 
Laboratory  and  c)  Ms.  Sigrid  Salo  ofNOAA. 


Sub-zero  temperatures,  strong  winds,  limited  daylight 
and  some  intimidating  wildlife  (Figure  2)  make  field¬ 
work  in  the  Arctic  a  difficult  task. 

A  majority  of  observations  of  the  Arctic  ice  cover, 
both  in  real-time  and  over  periods  of  decades,  have 
been  provided  by  satellite  imagery  and  drifting  buoys. 
Satellite  data  derived  from  the  Defense  Meteorological 
Satellite  Program  (DMSP)  Operational  Line-Scan  (OLS) 
Sensor,  the  National  Oceanic  and  Atmospheric 
Administration  (NOAA)  Advanced  Very  High 
Resolution  Radiometer  (AVHRR),  the  DMSP  Special 
Sensor  Microwave /Imager  (SSM/I),  and  the  Canadian 
RADARSAT  are  used  by  forecasting  centers  like  the 
U.S.  Navy/NOAA  National  Ice  Center  (NIC)  to  pro¬ 
vide  a  real-time  picture  of  conditions  in  the  Arctic.  A 
network  of  automatic  drifting  data  buoys  monitoring 
surface  air  temperature,  pressure  and  ice  motion,  has 
existed  in  the  Arctic  basin  since  early  1979,  first  under 
the  Arctic  Ocean  Buoy  Program  and  then  in  1991,  as 


part  of  the  International  Arctic  Buoy  Program  (IABP). 
These  data  have  been  used  to  support  real-time  opera¬ 
tions  in  the  Arctic  as  well  as  meteorological  and 
oceanographic  research  of  the  Arctic  basin.  More  infor¬ 
mation  on  the  IABP  is  available  at  http://iabp.apl. 
washington.edu  / . 

Forecasts  of  ice  conditions  are  produced  by  numer¬ 
ical  ice-ocean  models  that  use  these  observations  to 
help  specify  an  "initial"  state  and  then  run  forward  in 
time.  The  length  of  the  ice  forecasts  depends  on  the 
atmospheric  forcing  that  drives  the  model.  Usually,  the 
forcing  is  derived  from  an  atmospheric  forecast  model, 
and  extends  about  seven  days  into  the  future.  Longer 
forecasts  (to  30  days)  are  sometimes  generated  using 
persistent  atmospheric  conditions.  The  main  forecast 
products  that  describe  sea  ice  conditions  are  ice  thick¬ 
ness,  ice  concentration,  and  ice  drift. 

Real-time  sea  ice  observations,  analyses  and  fore¬ 
casts  are  now  available  from  ice  centers  around  the 
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Figure  2.  Polar  bears  inspecting  a  buoy  from  the  International  Arctic  Buoy  Program.  This  buoy  ceased  reporting  tempera¬ 
tures  shortly  after  this  picture  was  taken.  Photo  courtesy  ofD.G.  Barton,  US  Coast  Guard  (retired). 


world  (i.e.  Canadian,  Swedish,  Danish,  Japanese, 
Icelandic  centers).  These  centers  are  usually  responsi¬ 
ble  for  providing  information  on  ice  conditions  closest 
to  their  own  coastlines.  The  U.S.  NIC  is  actively  pro¬ 
viding  global  sea  ice  information  to  Navy  fleet  opera¬ 
tors  and  commercial  users.  In  addition  to  satellite 
observations,  aerial  ice  reconnaissance,  ship  /  shore  sta¬ 
tion  observations,  drifting  buoy  reports,  meteorologi¬ 
cal  guidance  products,  climatology  and  information 
from  international  partners,  the  NIC  uses  the  forecasts 
from  the  Navy's  sea  ice  forecasting  system,  the  Polar 
Ice  Prediction  System  (PIPS)  version  2.0  to  generate 
analyses  and  forecasts  of  sea  ice  conditions  in  the 
Arctic. 


The  Polar  Ice  Prediction  System  (PIPS) 

Since  the  late  1980s,  the  Naval  Research  Laboratory 
(NRL)  has  been  developing  sea  ice  forecast  systems  for 
operational  use  by  the  U.S.  Navy.  These  models  run 
operationally  on  computer  systems  at  the  Fleet 
Numerical  Meteorology  and  Oceanography  Center 
(FNMOC).  Products  are  disseminated  by  FNMOC  for 
use  at  the  NIC.  There  have  been  several  versions  of  the 
Polar  Ice  Prediction  System  in  the  past  that  have  led  up 
to  the  present  PIPS  2.0.  The  first  version  of  PIPS  (PIPS 
1.0),  operational  in  1987,  covered  the  Arctic  basin, 
Barents  and  Greenland  Seas  using  a  horizontal  grid  res¬ 
olution  of  127  km.  PIPS  1.0  was  followed  by  two  high¬ 
er  resolution  regional  models  for  the  Barents  Sea,  the 
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Figure  3.  The  PIPS  2.0  grid.  Every  fourth  grid  point  is 
plotted. 


Regional  Polar  Ice  Prediction  System — Barents  (RPIPS- 
B)  and  for  the  Greenland  Sea,  the  Regional  Polar  Ice 
Prediction  System — Greenland  Sea  (RPIPS-G).  The 
regional  models  differed  from  the  PIPS  1.0  model  in 
that  they  covered  a  smaller  area  and  used  higher  grid 
resolution  (25  km  vs.  127  km).  The  higher  resolution 
was  necessary  to  provide  a  better  estimate  of  ice  edge 
and  coastlines.  These  forecast  systems  made  use  of  an 
ocean  climatology  to  provide  ocean  stresses  and  heat 
fluxes  to  the  ice  model.  Atmospheric  stresses  and  flux¬ 
es  were  provided  to  the  ice  model  by  the  Navy 
Operational  Global  Atmospheric  Prediction  System 
(NOGAPS;  see  Rosmond  et  al.  in  this  issue).  Once  these 
forecast  systems  were  in  place,  additional  require¬ 
ments  were  identified.  High-resolution  forecasts  were 
needed  in  all  of  the  ice-covered  regions  of  the  Northern 
Hemisphere  and  sea  ice  models  required  coupling  to 
an  ocean  model  to  provide  daily  ocean  variability.  In 
response  to  these  requirements,  NRL  developed  the 
Polar  Ice  Prediction  System  2.0  (PIPS  2.0).  This  system 
became  operational  in  1996  and  replaced  PIPS  1.0, 
RPIPS-B  and  RPIPS-G. 

The  core  of  the  PIPS  2.0  is  a  coupled  ice-ocean 
model  that  consists  of  the  Hibler  ice  model  (Hibler, 
1979;  1980)  and  the  Bryan  and  Cox  ocean  model  (Cox, 
1984).  PIPS  2.0  coverage  extends  from  the  North  Pole 
south,  to  approximately  30°N  latitude  and  includes 
marginal  seas  such  as  the  Sea  of  Okhotsk,  the  Sea  of 
Japan  and  the  Yellow  Sea  in  the  Pacific  and  the  Gulf  of 
St.  Lawrence  and  the  Labrador  Sea  in  the  Atlantic 
(Cheng  and  Preller,  1996).  The  grid  resolution,  0.28°, 


varies  from  17-33  km  depending  on  the  location  of  the 
grid  square  within  the  spherical  coordinate  system 
(Figure  3).  PIPS  2.0  uses  a  rotated  coordinate  system  to 
avoid  the  problem  of  a  numerical  singularity  at  the  pole. 

For  both  the  ice  and  ocean  models,  the  lateral 
boundaries  are  defined  as  solid  walls  and  placed  away 
from  any  sea  ice  covered  regions  in  order  to  minimize 
their  effect  on  sea  ice  forecast  areas  of  interest. 
Following  the  technique  of  Sarmiento  and  Bryan 
(1982),  the  ocean  model  temperature  and  salinity  fields 
are  loosely  constrained  to  the  Levitus  (1982)  climato¬ 
logical  data  set  using  a  250  day  restoring  time  at  all  lev¬ 
els  of  the  ocean.  This  weak  constraint  on  the  ocean  tem¬ 
perature  and  salinity  does  not  adversely  affect  the  PIPS 
2.0  ability  to  forecast  ocean  variability  on  the  daily  to 
weekly  time  scales.  The  bathymetry  used  by  the  ocean 
model  is  derived  from  a  U.S.  Navy  5  minute  data  base 
called  the  Navy  Digital  Bathymetry  Data  Base  5'  x  5' 
(DBDBV;  NAVOCEANO,  1997). 

In  this  coupled  system,  the  ocean  model  provides 
input  to  the  ice  model  in  terms  of  ocean  currents,  salin¬ 
ity  and  heat  fluxes  (temperatures).  The  ice  model  pro¬ 
vides  salinity  and  temperature  changes  due  to  the 
growth  and  decay  of  sea  ice  and  surface  stresses  to  the 
ocean  model.  The  ocean  model  contains  15  vertical  lev¬ 
els,  each  level  increasing  in  thickness  with  depth.  Direct 
interaction  between  the  ice  and  ocean  model  occurs  in 
the  first  level  that  is  30  m  deep  and  represents  the  ocean 
mixed  layer. 

Figure  4  is  a  schematic  of  the  design  of  PIPS  2.0.  At 
the  center  of  the  system,  is  the  coupled  Hibler  ice 
model /Cox  ocean  model.  Atmospheric  forcing  is  pro¬ 
vided  by  NOGAPS.  The  atmospheric  forcing  fields 
used  by  the  forecast  system  are  surface  wind  stress,  sur¬ 
face  air  temperature,  surface  pressure,  surface  vapor 
pressure,  shortwave  radiation,  sensible  plus  latent  heat 
flux,  and  the  total  heat  flux.  The  system  produces  a  120- 
hour  forecast  each  day.  In  addition  to  the  output  of  the 
model  (ice  drift,  thickness,  concentration, 
growth /decay  of  ice,  ocean  currents,  temperature  and 
salinity),  restart  fields  consisting  of  the  model's  24-hour 
forecast  are  written  to  a  file  to  be  used  to  initialize  the 
next  day's  forecast.  If  the  model  restart  fields  are  not 
available,  a  model-derived  climatology  is  used  to  ini¬ 
tialize  the  forecast. 

The  quality  of  each  forecast  depends  strongly  on 
the  accuracy  of  the  initial  conditions  and  the  atmos¬ 
pheric  forecasts,  as  well  as  on  the  inherent  accuracy  of 
the  forecast  model.  If  the  unaltered  24-hour  forecast  is 
used  to  initialize  the  forecast,  these  errors  will 
inevitably  compound  eventually  leading  to  unaccept¬ 
ably  large  forecast  errors.  To  reduce  this  effect,  the  24- 
hour  forecast  fields  are  used  in  combination  with 
observed  ice  concentration  data  to  initialize  the  model. 
The  PIPS  2.0  assimilates  ice  concentration  data  from  the 
DMSP  Special  Sensor  Microwave /Imager  (SSM/I). 
These  data  are  chosen  for  three  reasons:  they  are  avail¬ 
able  daily,  in  real  time,  at  FNMOC;  the  resolution  of  the 
data  is  similar  to  the  resolution  of  the  model;  and  the 
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Polar  Ice  Prediction  System  2.0 


Figure  4.  A  schematic  of  the  PIPS  2.0  sea  ice  forecasting  system. 


data  cover  the  entire  model  domain.  The  SSM/ 1  bright¬ 
ness  temperatures  are  converted  into  ice  concentration 
data  using  the  Navy  CAL  / VAL  algorithm  (Hollinger  et 
al.,  1991)  at  FNMOC.  This  algorithm  is  sensitive  to  the 
ice / water  boundary  and  to  thin  ice.  As  such  it  provides 
a  good  estimate  of  the  location  of  the  ice  edge  and  a 
better  estimate  of  thin  ice  than  most  other  algorithms. 
However,  due  to  its  sensitivity  to  thin  ice,  it  often  satu¬ 
rates  too  quickly  to  100%  ice  concentration,  producing 
an  overestimate  of  high  ice  concentrations.  For  more 
information  on  algorithm  descriptions  and  compar¬ 
isons  see  http://www.natice.noaa. gov/science/prod¬ 
ucts  /descriptions  .html. 

After  the  PIPS  2.0  restart  fields  are  read  into  the 
model,  the  daily  SSM/I  ice  concentration  data,  inter¬ 
polated  to  the  model  grid,  are  read  in.  Based  on  the 
characteristics  of  the  CAL/VAL  algorithm,  the  SSM/I 
data  replaces  the  model-derived  ice  concentration  only 
at  locations  where  observed  concentration  is  greater 
than  80%  or  less  than  50%  and,  the  difference  between 
the  two  fields  is  greater  than  10%  or  5%  respectively. 
The  model-derived  ice  thickness  field  and  the  ocean 
surface  temperature  field  are  then  adjusted  to  be  con¬ 
sistent  with  the  concentration  data.  That  is,  if  the 


SSM/I  data  indicate  there  is  no  ice  where  the  model 
had  produced  ice,  ice  is  removed  from  the  model  field 
and  the  ocean  temperature  is  raised  one  degree  above 
freezing  to  restrict  immediate  ice  growth  in  this  loca¬ 
tion.  If  the  model  did  not  have  ice  in  a  location  where 
the  SSM/I  sees  ice,  a  small  amount  of  ice  (0.3-0. 6  m)  is 
added  to  the  model  thickness  field  and  the  ocean  tem¬ 
perature  is  set  to  the  freezing  temperature  of  sea  ice. 
These  adjustment  values  were  determined  from  a  series 
of  model  experiments.  As  the  adjustments  are  usually 
confined  to  a  few  grid  cells  near  the  ice  edge,  the  exper¬ 
iments  showed  no  serious  dynamical/numerical  prob¬ 
lems  associated  with  the  technique. 

Forecast  Products 

Since  the  end  of  July  1996,  PIPS  2.0  has  been  pro¬ 
ducing  operational  120-hour  forecasts  daily.  Each  day  a 
set  of  11  products  is  generated  from  the  PIPS  2.0  fore¬ 
cast  (Table  1).  These  products,  and  the  interval  at  which 
they  are  generated,  are  determined  by  the  NIC.  They 
are  transmitted  to  the  NIC  in  GRIB  format  and  may  be 
viewed  as  a  graphic  product  on  NIC  workstations.  The 
NIC  uses  PIPS  2.0  products  as  guidance  in  generating 
their  ice  analysis  and  forecast  products. 
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Table  1 

PIPS  2.0  Daily  Operational  Products 


Product 

Depth  *(m) 

Forecast  Time  (hrs) 

Ice  Thickness 

— 

0,  24,  48,  72,  96,  1 20 

Ice  Concentration 

— 

0,  24,  48,  72,  96,  1 20 

Ice  Drift 

— 

6,  12,  18,  24,  48,  72,  96,  120 

Divergence/ Convergence 

— 

6,  1 2,  1 8,  24,  48,  72,  96,  1 20 

Ocean  Currents 

15.0 

6,  12,  18,  24,  48,  72,  96,  120 

53.0 

6,  12,  18,  24,  48,  72,  96,  120 

110.0 

6,  12,  18,  24,  48,  72,  96,  120 

Ocean  Temperature 

15.0 

24,  48,  72,  96,  1 20 

53.0 

24,  48,  72,  96,  1 20 

110.0 

24,  48,  72,  96,  1 20 

Ocean  Salinity 

15.0 

24,  48,  72,  96,  1 20 

‘Depth  is  at  the  center  of  each  level 


Figures  5b-d  are  examples  of  standard  products 
from  the  system:  ice  displacement  (based  on  the  ice 
drift),  ice  concentration  and  ice  thickness.  Figure  5a 
represents  the  observed  ice  motion  derived  from  the 
available  IABP  drifting  buoys.  Figures  5a  and  5b  show 
the  qualitative  similarities  between  the  PIPS  2.0  fore¬ 
cast  and  the  observed  ice  motion.  PIPS  2.0  also  has  the 
capability  to  forecast  ocean  currents,  ocean  tempera¬ 
ture  and  salinity.  Recently  these  fields  were  requested 
by  the  NIC  as  input  to  a  marginal  ice  zone  specific  ice 
model  that  is  being  developed  and  tested  for  the 
center's  use. 

Future  Sea  Ice  Forecast  Systems— PIPS  3.0 

In  the  pre-operational  evaluation  of  PIPS  2.0  prod¬ 
ucts  (Preller  and  Posey,  1995)  the  original  PIPS  and 
PIPS  2.0  ice  drift  were  evaluated  against  Arctic  buoy 
data.  These  results  showed  that  PIPS  2.0  produced  an 
improved  ice  drift  forecast  over  PIPS  with  a  decrease  in 
the  Root  Mean  Squared  (RMS)  Error  of  about  10-15%. 
Due  to  the  combination  of  the  increased  resolution  and 
the  incorporation  of  an  ocean  model  versus  an  ocean 
climatology,  the  PIPS  2.0  ice  edge  was  a  substantial 
improvement  over  the  existing  PIPS  models,  including 
the  higher  resolution  regional  models. 

A  recent  study  by  a  group  of  scientists  from  the 
NIC  and  NOAA  (Van  Woert  et  al.,  2001),  showed  that 
although  the  PIPS  2.0  forecasts  (48-hour)  were  better 
than  persistence  on  average,  there  were  still  substantial 
biases  in  its  prediction  of  the  growth  and  decay  of  sea 
ice  in  the  marginal  ice  zone.  PIPS  2.0  often  over-pre- 
dicts  the  amount  of  ice  in  the  Barents  Sea  and  therefore 
often  places  the  ice  edge  too  far  south.  In  contrast,  PIPS 
2.0  often  under-predicts  the  ice  extent  in  the  Labrador 
Sea  and  Hudson  Bay.  As  forecasts  are  a  combination  of 


many  different  components  (Figure  4) — atmospheric 
forcing,  oceanographic  forcing,  model  parameteriza- 
tions  and  initialization — the  combination  of  even  small 
errors  from  each  of  these  pieces  could  add  up  to  a  sub¬ 
stantial  error  in  the  forecast. 

In  recent  years,  there  has  been  a  request  for  the 
capability  to  predict  the  location  of  regions  of  open 
water  within  the  ice  pack.  Approximately  1%  of  the 
central  Arctic  ice  pack  is  estimated  to  be  open  water  in 
winter.  In  summer,  the  amount  increases  to  10-20% 
(Gow  and  Tucker,  1990).  Openings  in  the  ice  pack  usu¬ 
ally  appear  in  two  forms:  leads  and  polynyas.  Leads  are 
crack-like  openings  in  the  ice  pack,  usually  tens  to  hun¬ 
dreds  of  meters  wide,  while  polynyas  are  larger,  more 
permanent  openings.  Leads  form  in  regions  where 
large-scale  divergent  wind  patterns  produce  divergent 
stresses  in  the  ice  causing  the  ice  to  break  apart.  Leads 
appear  as  linear  features  in  compact  ice  that  often  have 
a  preferential  orientation  (direction).  Polynyas  form 
near  coastlines  where  divergence  or  shearing  of  the 
pack  ice  separates  the  ice  from  the  coast  or  from  a  fast 
ice  boundary.  Leads  and  polynyas  are  significant 
because  they  provide  traversable  areas  in  the  ice  pack 
for  ships.  In  addition  they  are  important  because  they 
are  regions  of  very  large  heat-loss  to  the  atmosphere. 
Winter  leads  rapidly  refreeze,  due  to  the  heat  loss  to  the 
atmosphere,  forming  young  thin  ice  that  is  easily 
deformed.  Leads  and  polynyas  are  a  major  source  of 
new  ice  growth  in  winter.  Although  PIPS  2.0  can  pre¬ 
dict  large-scale  polynyas,  it  does  not  have  the  capabili¬ 
ty  to  produce  smaller  polynyas  or  leads  or  to  provide 
guidance  on  lead  orientation. 

PIPS  2.0  is  based  on  the  scientific  research  efforts  of 
late  1970s  and  early  1980s.  During  the  past  10  years 
(1990s),  great  strides  have  been  made  in  understanding 
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sea  ice  dynamics  and  thermodynamics  as  well  as 
observing  ice  conditions.  An  additional  and  very 
important  factor  in  the  improvement  of  ice  modeling 
and  forecast  capabilities  is  the  advance  in  computer 
technology  over  the  past  10  years.  Computer  codes 
now  make  use  of  multiple  processors  and  can  perform 
more  extensive  computations  in  operationally  accept¬ 
able  time  periods. 

In  1998,  the  Office  of  Naval  Research  (ONR)  and 
the  Oceanographer  of  the  Navy  via  the  Space  and 


Naval  Warfare  Systems  Command  (SPAWAR)  joined 
forces  to  fund  an  effort  to  combine  this  new  technology 
and  data  into  an  improved  sea  ice  forecasting  system. 
This  system,  aptly  named  PIPS  3.0  (http://www.oc.- 
nps.navy.mil/~pips3)  is  presently  being  developed 
through  a  joint  effort  among  the  U.S.  Naval 
Postgraduate  School,  various  other  academic  institu¬ 
tions  and  the  Naval  Research  Laboratory.  The  PIPS  3.0 
ice  model  will  contain  an  improved  calculation  of  ice 
growth/ decay  based  on  the  use  of  a  multi-level  ice 


Figure  5.  IABP  buoy  drift  and  PIPS  2.0  72-hour  forecast  products  from  December  30,  2000.  a)  IABP  buoy  drift  in  cm/sec  and 
sea  level  pressure  in  millibars,  b)  72-hour  ice  displacement  in  kilometers,  c)  ice  concentration  in  percent  and  d)  ice  thickness 
in  meters. 
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Figure  6.  Zn  degree,  high  resolution,  con-pled  model  domain  and  bathymetry.  Contours  of  25,  50, 100,  500, 1000,  2000,  3000, 
4000,  5000  and  6000  m  isobaths  are  shown.  More  information  about  the  model  is  available  at 

http://www.ocmps.navy.mil/~pips3. 


thickness  formulation.  Ice  edge  forecasts  will  be 
improved  by  using  higher  horizontal  resolution.  Ice 
motion  and  ice  edge  location  will  be  improved  by  the 
assimilation  of  satellite  derived  ice  drift  data.  In  addi¬ 
tion,  higher  resolution  and  more  realistic  ice  rheology 
will  improve  the  PIPS  ability  to  predict  areas  of  lead 
formation  and  lead  orientation. 

These  improvements  to  the  PIPS  2.0  are  being  test¬ 
ed  in  an  incremental  fashion.  The  first  tests  were 
designed  to  examine  the  effects  of  higher  resolution  on 
an  ice  model  similar  to  that  of  PIPS  2.0  coupled  to  an 
improved  ocean  model  with  improved  bathymetry.  At 
the  same  time,  techniques  for  the  assimilation  of  satel¬ 
lite  derived  ice  drift  into  the  ice  model  were  developed 
and  tested.  Results  from  these  tests  are  described 
below. 

High  Resolution  Coupled  Model  Tests 

Similar  to  the  PIPS  2.0,  the  coupled  ice-ocean 
model  used  in  these  experiments  extends  from  the 
North  Pacific  (at  ~30°N)  through  the  Arctic  Ocean  and 
the  Nordic  Seas  into  the  North  Atlantic  (to  ~45°N; 
Figure  6)  and  uses  a  rotated  coordinate  system.  The  ice 
and  ocean  model  horizontal  grid  is  configured  at  Zii 
(~9  km)  with  45  vertical  levels  in  the  ocean.  In  the  ver¬ 
tical  direction,  the  upper  100  meters  of  the  ocean  are 


resolved  by  eleven  layers,  and  the  upper  500  meters  of 
the  ocean  by  nineteen  layers.  This  resolution  allows 
detailed  representation  of  bathymetry  over  the  vast 
Arctic  shelves  and  in  marginal  seas  as  well  as  many 
local  bathymetric  features  important  to  large-scale 
ocean  circulation  (e.g.  St.  Anna  Trough,  Canadian 
Archipelago,  or  Bering  Strait).  The  bathymetry  data 
used  to  define  the  grid  shown  in  Figure  6  consists  of  the 
ETOP05  database  (National  Geophysical  Data  Center, 
1988),  NRL  charts  and  Canadian  Hydrographic  Center 
charts  for  latitudes  south  of  64°N.  North  of  64°N,  the 
2.5  km  resolution  International  Bathymetric  Chart  of 
the  Arctic  Ocean  (IBCAO)  database  (Jakobsson  et  al., 
2000)  is  used.  These  data  represent  a  considerable 
improvement  over  previous  data  available  for  this  area 
used  in  PIPS  2.0.  This  high  resolution,  coupled  ice- 
ocean  model  is  an  extension  of  previous  modeling 
research  using  a  horizontal  resolution  of  Zb  and  30  ver¬ 
tical  levels  (Maslowski  et  al.,  2000,  2002;  Zhang  et  al., 
1999,  http://www.oc.nps.navy.mil/sbi).  The  coupled 
model  adapts  the  Los  Alamos  National  Laboratory 
(LANL)  global  Parallel  Ocean  Program  (POP)  ocean 
model  with  a  free  surface  (Maltrud  et  al.,  1998; 
Dukowicz  and  Smith,  1994).  The  sea  ice  model  at  pres¬ 
ent  is  similar  to  PIPS  2.0  with  a  viscous-plastic  rheolo¬ 
gy,  the  zero-layer  approximation  of  heat  conduction 
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through  ice  and  a  surface  energy  budget  calculated  fol¬ 
lowing  Parkinson  and  Washington  (1979).  The  main 
difference  between  this  ice  model  and  the  PIPS  2.0  ice 
model  is  its  ability  to  utilize  modern  parallel  computer 
architectures. 

A  27-year  spin-up  integration  of  the  coupled  ice- 
ocean  model  has  been  run,  forced  with  a  daily-aver¬ 
aged  annual  cycle  of  climatological  atmospheric  fields 
derived  from  the  European  Centre  for  Medium-range 
Weather  Forecasts  (ECMWF)  1979-1993  reanalysis. 
Forcing  fields  include:  10-meter  winds,  surface  pres¬ 
sure,  surface  air  temperature  and  dew  point,  and  long¬ 
wave  and  short-wave  radiation.  During  the  spinup,  the 
ocean  surface  temperature  and  salinity  are  restored 
with  a  30-day  time  scale  to  monthly  mean  climatology 
derived  from  the  University  of  Washington's  Polar 
Science  Center  Hydrographic  Climatology  (PHC; 
Steele  et  al.,  2001). 

Following  the  27-year  integration,  an  additional 
12-year  run  has  been  completed  using  the  repeated 
1979  ECMWF  annual  cycle  for  the  first  6  years  and 
1979-1981  inter-annual  fields  for  the  last  6  years.  This 
approach  has  been  used  to  force  the  sea  ice  and  ocean 
states  towards  conditions  of  the  late  1970s  and  early 
1980s.  These  conditions  may  then  be  used  to  initialize 


the  new  version  of  the  PIPS  2.0  ice  model  and  spin  it  up 
to  the  present  day  for  forecast  use. 

Understanding  the  ice  thickness  distribution  and 
its  variability  is  important  in  both  short-term  opera¬ 
tional  forecasts  and  longer-term  climate  studies. 
Model-derived  ice  thickness  fields  are  useful  in  obtain¬ 
ing  a  large-scale,  but  detailed  picture  of  ice  distribution 
and,  more  importantly,  an  understanding  of  its  tempo¬ 
ral  variability. 

In  Figure  7,  the  3-year  mean  ice  thickness  distribu¬ 
tion  is  presented  from  the  first  three  years  of  the 
1979-1981  output.  The  modeled  distribution  of  ice 
thickness  compares  reasonably  well  to  that  known 
from  observations  (Bourke  and  Garrett,  1987;  Bourke 
and  McClaren,  1992).  In  agreement  with  data,  the  thick¬ 
est  ice  (>  6.0  m)  in  the  Arctic  is  found  along  the 
Canadian  Archipelago  and  northern  coast  of 
Greenland.  Farther  north  near  the  pole,  the  ice  thick¬ 
ness  decreases  to  values  of  3. 0-3. 5  m.  This  simulation 
also  shows  relatively  thick  ice  (>  3.5  m)  on  the  East 
Siberian  shelf.  Since  few  observations  are  available  for 
that  region  and  that  time,  it  is  difficult  to  quantitatively 
verify  these  results.  In  the  1990s,  measurements  of  ice 
thickness  were  made  in  the  central  Arctic  Ocean  from 
U.S.  Navy  submarines  and  at  a  number  of  key  point 
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locations  (e.g.  the  Fram  Strait)  using  upward  looking 
sonar.  Model  results  will  be  verified  against  these 
recent  observations  once  the  model  output  for  that 
time  period  becomes  available.  Based  on  the  compar¬ 
isons  made  thus  far,  these  model  results  look  promis¬ 
ing  in  terms  of  forecasting  ice  thickness  and  ice  edge 
variability. 

Coupled  ice-ocean  models,  such  as  those  used  by 
PIPS  2.0  and  PIPS  3.0,  are  most  often  applied  in  aca¬ 
demic  studies  of  global  climate  issues  such  as  climate 
change.  In  these  studies  the  models  are  used  for 
decadal  or  longer  simulations  and  focus  on  long  term 
changes  in  the  overall  pattern  of  ice  thickness  and  ice 
extent.  This  application  of  the  coupled  model  requires 
access  to  powerful  computers  and  large  quantities  of 
computer  time.  However,  coupled  ice-ocean  models 
used  for  forecasting  face  a  different  set  of  require¬ 
ments.  As  forecasts  are  concerned  with  sea  ice  variabil¬ 
ity  on  much  shorter  time  scales,  the  models  must  be 
designed  to  produce  the  most  accurate  daily  changes 
in  ice  concentration,  ice  edge  location  and  ice  motion. 
These  models  are  required  to  go  through  rigorous 
validation  studies  to  prove  their  capability  to  produce 
accurate  short  term  variability.  Data  assimilation  plays 
a  major  role  in  the  accuracy  of  these  forecasts.  Once 
operational,  continuous  quality  control  and  evaluation 
of  the  products  may  be  used  to  upgrade  the  system  and 
improve  forecast  accuracy.  In  addition,  forecast 
systems  are  limited  by  the  amount  of  computer 
resources  available  for  each  forecast  as  they  compete 
with  other  forecast  models  each  day.  Therefore  the 
combined  model/assimilation  system  must  be 
designed  to  "fit"  within  these  limitations.  This  places 
restrictions  on  the  "size"  or  grid  resolution  of  the  mod¬ 
els  as  well  as  the  complexity  of  the  model  parameteri- 
zations  and  the  data  assimilation  techniques.  Each  of 
these  issues  must  be  taken  into  account  when  develop¬ 
ing  a  new  forecast  system. 

PIPS  3.0  Ice  Motion  Data  Assimilation 

As  stated  earlier,  a  popular  strategy  to  reduce  the 
initial  error  in  numerical  forecast  models  is  to  use 
information  from  real-time  observations  to  correct  the 
previous  forecast  field  through  data  assimilation.  Data 
assimilation  analysis  makes  use  of  a  predictor-correc¬ 
tor  strategy,  where  the  model  predicted  ice  state 
(known  as  the  background  field)  is  corrected  by  what¬ 
ever  observations  are  available.  The  specific  process  of 
correction  fundamentally  determines  the  nature  of  the 
assimilation.  Current  ice  state  observations  tend  to 
be  either  scattered  through  space  or  infrequent  in  time. 
They  may  also  contain  noise  that  can  introduce 
non-physical  structure  into  the  model.  The  direct 
replacement  of  model  values  with  the  observed  values 
(known  as  insertion)  can  have  undesirable  conse¬ 
quences.  It  is  therefore,  preferable  to  employ  some 
sort  of  self-consistent  procedure  that  reduces  the  error 
and  distributes  the  observational  knowledge  widely 


over  the  domain.  Statistical  assimilation  methods  meet 
these  criteria. 

One  of  the  first  investigations  of  the  assimilation  of 
ice  motion  data  into  a  stand-alone  ice  model  (similar  to 
the  original  PIPS  1.0)  employs  the  optimal  interpolation 
methodology  (Meier  et  al.,  2000).  Optimal  interpolation 
is  the  simplest  type  of  statistical  assimilation.  Statistical 
assimilation  means  that  the  background  field  is  correct¬ 
ed  to  "minimize  the  error  statistics"  of  the  assimilated 
field.  This  can  be  done,  only  if  the  error  statistics  of  the 
observations  and  the  model  are  known.  The  technique 
is  "optimal"  since  the  error  variance  is  minimized. 

Individual  ice  velocity  vectors  of  the  background 
field  are  corrected  by  applying  a  linear  combination  of 
the  nearby  observed  ice  motions.  The  correction  coeffi¬ 
cient  for  each  observed  value  results  from  solving  a  lin¬ 
ear  system  dependent,  among  other  things,  on  the  ratio 
of  the  model  to  data  error  covariances.  When  these 
weights  are  applied  to  the  corresponding  observations, 
the  resulting  correction  minimizes  the  error  variance  of 
the  corrected  solution.  This  formula  illustrates  that  it  is 
not  sufficient  to  be  concerned  with  the  observed  ice 
motions  alone.  The  error  statistics  for  the  model  and  the 
data  are  just  as  important.  They  instruct  the  algorithm 
as  to  which  observations  to  use,  and  what  relative 
weight  they  must  have,  in  order  for  the  correction  to 
minimize  the  error. 

Two  types  of  ice  motion  observations  are  used.  The 
first,  which  provides  the  data  to  be  assimilated,  is  not 
directly  observable.  It  is  derived  from  the  SSM/I  pas¬ 
sive  microwave  brightness  temperature  imagery 
through  a  technique  known  as  feature  tracking. 
Individual  swaths  of  the  SSM/I  satellite  are  combined 
or  "composited"  into  a  daily  gridded  field.  The  dis¬ 
placement  of  features  common  to  the  fields  on  consec¬ 
utive  days  can  provide  an  estimate  of  the  ice  velocity. 
The  derived  ice  motion  observations  used  here  were 
obtained  from  the  Polar  Remote  Sensing  Group  at  the 
NASA  Jet  Propulsion  Laboratory  (Kwok  et  al.,  1998). 
They  were  chosen  for  the  assimilation  analysis  because 
they  are  similar  to  the  real-time  ice  motion  products 
that  will  be  available  from  FNMOC.  Ice  motion  obser¬ 
vations  are  available  daily  from  the  beginning  of 
October  through  the  end  of  May.  The  data  are  not  avail¬ 
able  for  the  summer  months  due  to  the  susceptibility  of 
the  passive  microwave  signal  to  excessive  error  from 
moisture  in  the  form  of  either  heavy  summer  cloud 
cover  or  the  formation  of  melt  ponds  in  the  ice. 

The  second  set  of  observations  comes  from  IABP 
drifting  ice  buoys  (Rigor  and  Heiberg,  1997).  In  the 
operational  configuration  of  the  model,  the  buoy 
motions  will  be  included  as  observations  to  be  assimi¬ 
lated.  This  will  be  especially  important  in  the  summer 
months  when  SSM/I  derived  motions  are  not  available. 
The  buoys  were  withheld  in  this  study  to  be  used  for 
validation  purposes,  and  to  compute  the  error  covari¬ 
ances.  For  the  period  of  the  study,  there  are  roughly  30 
buoys,  irregularly  distributed  throughout  the  Arctic,  at 
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Figure  8.  Ice  drift  for  the  IABP  drifting  buoys  (red),  the  model  (green),  and  the  model  with  assimilation  (black). 


any  time.  Each  buoy  transmits  its  location  twice  per 
day  The  ice  buoy  velocity  can  then  be  calculated  from 
the  spatial  displacement  of  the  buoy  location. 

As  a  test  of  concept,  this  proposed  PIPS  3.0  ice 
motion  assimilation  method  was  instituted  in  a  previ¬ 
ously  tested  high  resolution  fully  coupled  sea- 
ice/ocean  model  (Zhang  et  al.,  1999).  The  coupled  sea- 
ice/ocean  model  employs  a  reduced  domain 
consisting  of  the  Arctic  basin,  with  a  closed  Bering 
Strait,  and  extending  south  into  the  Atlantic  to  about 
50°N.  The  spatial  resolution  of  the  model  is  roughly  18 
km:  half  the  resolution  proposed  for  the  PIPS  3.0 
model.  The  model  can  be  run  in  either  control  (no 
assimilation)  or  assimilation  mode. 

Two  model  runs  were  generated  for  the  period 
from  January  1992  through  December  1994.  The 
ECMWF  reanalysis  product  provided  the  atmospheric 


forcing  in  both  configurations.  In  the  control  run,  no 
assimilation  takes  place.  In  the  assimilative  run,  the 
model  digests  the  derived  SSM/I  ice  motion  product 
for  the  months  of  October  through  May  of  each  year. 
During  the  summer  months,  no  assimilation  takes 
place.  The  drifting  ice  buoy  data  are  withheld  as  vali¬ 
dation.  Figure  8  shows  selected  trajectories  of  the  IABP 
drifting  buoys,  the  control  run  and  the  assimilative  run. 
Trajectories  from  the  control  run  and  the  assimilative 
run  are  produced  by  integrating  forward  in  time  from 
the  initial  buoy  location  using  a  standard  fourth-order 
Runge-Kutta  method  and  the  respective  daily  instanta¬ 
neous  velocity  field.  The  trajectories  displayed  repre¬ 
sent  the  range  of  ice  motion  behavior  seen  during  the 
spring  of  1992,  and  the  period  covering  the  fall  of  1992 
through  the  spring  of  1993.  The  subset  of  trajectories 
displayed  in  the  figure  was  chosen  for  its  coverage  of 
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the  figure  domain.  The  trajectories  for  the  control  run 
and  the  assimilative  run  have  the  same  starting  loca¬ 
tions,  and  the  same  starting  and  ending  dates  as  the 
buoys. 

Consider  the  trajectories  numbered  1-10.  Both  the 
control  run  (green)  and  the  assimilative  run  (black) 
generated  motions  are  consistently  slower  than  the 
buoy  motion  (red).  The  pattern  of  buoy  twists  and 
turns  is  strongly  reflected  in  the  assimilative  run  trajec¬ 
tories  for  about  a  half  of  the  cases.  The  error,  for  the 
assimilative  run  trajectories,  is  largely  due  to  the  dif¬ 
ference  in  speed.  The  error  in  the  control  run  is  mani¬ 
fested  in  terms  of  an  even  greater  speed  differential  as 
well  as  large  errors  in  the  direction  of  the  motion.  In 
nine  out  of  the  ten  cases,  the  assimilative  run  produces 
a  trajectory  that  more  closely  follows  the  buoy  trajecto¬ 
ry  then  does  the  control  run's  trajectory.  The  assimila¬ 
tion  process  generally  increases  the  ice  speed,  and  pre¬ 
serves  the  observed  direction  of  the  ice  motion.  On 
average,  the  trajectories  with  assimilation  contain  half 
the  error  of  the  trajectories  without  assimilation. 

The  remaining  two  trajectories  (11  and  12)  illus¬ 
trate  a  different  situation.  The  marginal  ice  zone  (MIZ) 
in  the  Greenland  Sea  is  extremely  dynamic,  commonly 
exhibiting  the  largest  ice  speeds  of  anywhere  in  the 
domain.  In  addition,  south  of  80°N  there  are  very  few 
SSM/I  derived  observations  that  can  be  assimilated, 
and  those  that  exist  tend  to  have  large  error  bounds.  As 
a  result,  trajectories  south  of  80°N  can  have  quite  large 
errors.  For  trajectory  12,  which  starts  east  of 
Spitsbergen,  the  buoy  moves  clockwise,  encircling  the 
island.  In  contrast,  the  trajectories  representing  the 
control  run  and  the  assimilative  run  speed  south  into 
the  Norwegian  Sea.  The  similarity  between  these  two 
trajectories  is  due  to  the  lack  of  SSM/I  derived  ice 
motions  to  assimilate.  The  difference  in  the  path  of  the 
two  trajectories  reflects  the  difference  between  the  ice 
states  for  that  location  in  the  two  model  runs.  Since  the 
current  ice  state  for  a  given  location  is  a  function  of  its 
history,  each  model  can  have  ice  with  different  concen¬ 
tration  and  thickness  characteristics.  The  different  ice 
characteristics  can  result  in  different  reactions  to  the 
same  forcing. 

Trajectory  11  shows  the  buoy  flowing  south  into 
the  Fram  Strait  and  along  eastern  Greenland.  Initially, 
the  three  trajectories  are  indistinguishable.  As  they 
approach  Fram  Strait,  the  trajectory  from  the  control 
run  begins  to  diverge  from  the  path  of  the  other  two. 
From  the  starting  location  to  approximately  79°N,  the 
assimilative  run  produces  a  trajectory  that  coincides 
well  with  the  buoy.  Past  this  point,  where  there  are  no 
longer  any  SSM/I  derived  motions  for  the  model  to 
assimilate,  the  three  trajectories  rapidly  diverge. 
Eventually,  both  the  control  run  and  the  assimilative 
run  produce  trajectories  that  travel  up  to  twice  as  far  as 
the  buoy,  with  the  control  run's  trajectory  traveling  the 
farthest.  Neither  of  the  model  runs,  with  or  without 
data  assimilation,  realistically  represents  the  buoy 
motion  south  of  the  Fram  Strait.  These  two  cases  illus¬ 


trate  the  need  for  the  assimilation  of  the  buoy  motions 
in  regions  where  the  SSM/I  derived  motions  are 
not  available.  In  addition,  the  buoy  tracks  can  provide 
badly  needed  data  during  the  summer  season, 
when  SSM/I  derived  products  are  not  available. 
Further  details  on  the  assimilative  model  can  be  found 
in  Stark  (2001). 

Summary 

The  U.S.  Navy  has  been  running  sea  ice  forecast 
systems  operationally  since  the  late  1980s.  These  mod¬ 
els  have  been  gradually  upgraded  through  the  years  to 
higher  resolution,  data  assimilative  models  that  span 
more  of  the  sea  ice  covered  oceans  in  the  Northern 
Hemisphere. 

The  models  have  provided  forecasts  of  ice  motion, 
ice  thickness  and  ice  coverage  over  regions  that  are 
physically  difficult  to  observe  due  to  hostile  environ¬ 
mental  conditions.  Advancements  in  satellite  technolo¬ 
gy  provide  these  forecast  systems  with  continuously 
improving  data  for  assimilation,  thus  enabling 
improvements  to  the  daily  forecasts.  Advancements 
in  computer  technology  provide  these  forecast 
systems  larger  and  faster  computers  to  generate 
improved  forecasts. 

The  most  recent  upgrade  to  the  Navy's  ice  predic¬ 
tion  capability  is  the  development  of  the  next  genera¬ 
tion  forecast  system,  PIPS  3.0.  Improvements  to  this 
new  forecast  system  include,  higher  horizontal  resolu¬ 
tion,  a  more  sophisticated  ocean  model,  improved  data 
assimilation  and  perhaps  most  important,  an  improved 
sea  ice  model.  The  ongoing  upgrade  of  the  sea  ice 
model  will  include  a  Lagrangian  formulation  for  calcu¬ 
lating  a  multi-category  ice  thickness  distribution,  a 
snow  layer,  a  brine  pocket  parameterization,  non-linear 
profiles  of  temperature  and  salinity  (Bitz  and 
Lipscomb,  1999),  and  a  Coulombic  yield  curve  for  the 
viscous-plastic  rheology  (Hibler,  2000).  These  improve¬ 
ments  are  geared  towards  providing  better  forecasts  of 
ice  edge,  ice  motion,  ice  thickness  and  regions  of  lead 
formation  and  lead  orientation.  The  PIPS  3.0  is  present¬ 
ly  going  through  its  final  development  and  will  begin 
its  adaptation  for  operational  use  next  year  with  a 
scheduled  transition  into  operational  use  in  late  2002  or 
early  2003.  Ei2 
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